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A  RANKINE-HUGONIOT  EMULATING  (DENSITY  vs.  VELOCITY) 
RELATIONSHIP  FOR  CFD  USAGE  WITHIN  AN  INVISCID  SHOCK  WAVE 


by 

B.W.B.  SHAW 
SUMMARY 

To  enable  the  Full  Potential  Equation  method  of  computing  transonic  inviscid  flows 
to  accurately  predict  the  conditions  at  exit  from  a  shock  wave  -  a  desirable  aim  which  the 
method  fails  to  achieve  to  a  greater  or  lesser  extent  -  an  improved  relationship  has  been 
derived  between  the  density  (p)  and  velocity  (q)  inside  an  inviscid  shock  wave.  This 
relationship  replaces  the  conventional  isentropic,  isenergic  relationship  normally  applied 
there.  The  shock  exit  conditions  thereby  obtained  are  the  correct  Rankine-Hugoniot  values, 
for  all  shock  inlet  Mach  numbers. 

Zonal  application  of  this  improved  (p,q)  shock  wave  internal  relationship  should 
eliminate  the  prime  cause  of  solution  inaccuracy  when  using  the  conservative  Full  Potential 
Equation  method. 

Two  alternative  versions  of  this  relationship  are  provided;  CFD  usage  will  indicate 
which  is  superior. 

The  analysis  has  been  applied  first  to  normal  shock  waves  and  then  extended  to 
oblique  shock  waves. 
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1  INTRODUCTION 


The  abrupt  changes  occurring  across  a  physical  shock  wave  are  not  truly  isen- 
tropic.  Instead,  they  obey  the  Rankine-Hugoniot  relationship  and  depart  more 
and  more  from  being  isentropic  flow  changes  as  the  shock  wave  strength  in¬ 
creases.  that  is  as  the  shock  wave  entry  Mach  number  Min  increases.  How¬ 
ever  the  full  potential  equation  (FPE)  approach,  which  has  been  widely  used 
(e.g.  Refs.  1-8)  in  CFD  solutions  of  the  transonic  flows  past  aerofoils,  wings,  bod¬ 
ies.  wing-body  combinations,  intakes,  cascades  and  helicopter  blades,  makes  the 
approximation  that  changes  of  state  are  iseuiropic  throughout  the  entire  flow 
field  including  any  shock  wave(s)  present.  Thus  when  used  with  shocked  flows 
the  FPE  flow  model  may  well  introduce  significant  errors  (see  Refs.  9.  10)  in  the 
computed  physical  quantities  of  the  flow  and  thereby  in  the  global  quantities 
Cl,  Cd  etc.  For  air  flowing  through  a  normal  shock  wave  these  physical  errors 
are  listed  in  Table  1  and  illustrated  at  Figs.  2-4.  (In  Table  1  AM<.z  is  the 
error  in  shock  wave  exit  Mach  number  Me«,  APR  and  ADR  are  the  percentage 
errors  in  the  pressure  ratio  PR  and  density  ratio  DR  across  the  shock  wave.) 

These  Table  1  errors  vary  approximately  as  a  power  of  (Min^  -  1).  and  as 


Min 

apr9^o 

ADR% 

1.0 

0 

0 

0 

1.1 

-.OOr 

1 

i 

1.3 

-.043 

6i 

5 

1.5 

-.09, 

16 

13i 

1.6 

-.Us 

22i 

19 

Table  1:  ERRORS  RESULTING  FROM  THE  FPE  APPROXIMATION  OF  A 
NORMAL  SHOCK  WAVE  (7=  1-4) 


Mjn  increases  from  1.0  they  are  at  first  insignificant  but  in  due  course  increase 
rapidly  to  large  values.  There  is  thus  a  need  for  an  improved  model  of  the  flow 
changes  across  a  shock  wave,  for  use  with  the  FPE  approach,  particularly  for 
transonic  flows  with  relatively  high  Mach  numbers  at  the  shock  wave  entry. 

The  mathematical  formulation  of  the  conservative  FPE  flow  model  consists  of 
a  set  of  equations  in  the  FPE  dependent  variables  —  potential  function  #. 
velocity  q  (together  with  its  components)  and  density  p.  In  particular,  it  in¬ 
volves  the  use  of  an  algebraic  (p  vs.  q)  equation  that  is  based,  as  previously 
indicated,  on  the  assumption  of  universally  isentropic  flow;  this  is  the  source  of 
the  errors  mentioned  above.  The  analysis  which  follows  describes  the  develop¬ 
ment  of  an  alternative  (p  vs.  q)  continuous  relationship  that  is  non-isentropic 
and  that  ensures  that  shock  wave  exit  conditions  are  as  given  by  the  Rankine- 
Hugoniot  relations.  The  relationship  applies  only  in  the  entire  interior  region  of 
an  "invtscid^  shock  wave  wafer.  *  It  is  intended  for  zonal  use  with  FPE  flow 

'The  phrase  -shock  wave  wafer"  is  used  here  and  subsequently  to  emphasiie  the  small  but 
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comput.itions  and  should  eliminate  the  errors  referred  to  above.  Elsewhere 
in  the  flow  the  usual  FPE  isentropic.  isenergic  flow  law  giving  {f>  vs.q)  woirld 
apply,  though  with  an  allowance  downstream  of  the  shock  wave  for  the  reduced 
stagnation  density  there. 

In  fashioning  this  relationship  the  principal  aims  have  been  that: 

•  The  relationship  should  ensure  the  correct  —  i.e.  Rankine-Hugoniot  — 
values  of  velocity  q  and  density  p  at  e.'dt  from  the  shock  wave  wafer. 
Consequently  mass,  energy  and  total  impulse  must  be  conserved  in  the 
transition  from  entry  to  exit  locations  of  the  shock  wave  wafer. 

•  The  flow  quantities  must  be  piecewise  continuous  across  the  shock  wave 
wafer  entry  face. 

It  also  appeared  desirable  that: 

•  p  should  be  a  monotonic  function  of  q  (to  accord  with  physical  reality), 
and 

•  The  finally  chosen  algebraic  expression  for  p  =  p(q)  should  involve  a 
minimum  of  computational  effort,  consistent  with  satisfying  the  first  two 
requirements  with  high  accuracy, 

7  has  been  assumed  to  be  1,4. 


2  CONDITIONS  AT  EXIT  FROM  A  NORMAL 
SHOCK  WAVE 


The  procedure  adopted  has  been  to  derive  first  a  (p,q)  relationship  for  a  normal 
shock  wave.  In  the  .Appendix  this  has  then  been  extended  to  the  more  general 
case  of  an  oblique  shock  wave. 

Using  q.  p,  p  and  a  to  denote  speed,  density,  pressure  and  speed  of  sound, 
the  continuity,  momentum  and  energy-conservation  equations  across  a  normal 
shock  wave,  with  upstream  and  downstream  conditions  denoted  by  suffices  "I" 
and  “2”  (Fig.  1  refers),  are  as  foliows:- 


continuity 


P2^2  =  Pl<7l 


(1) 


momentum 


P29I -I-P2  =  P\q\  +  Pi 


i.e. 

P29^[1  4-  alK'tql)]  =  p,qf[l  -t-  a\l('iq\)\ 


finite  thickness  of  computational  (as  well  as  physical)  shock  waves. 


'2' 


'1' 


{Qi'  Pl*Pl*®l) 


Figure  1:  NORMAL  SHOCKWAVE  STATES 


energy 

i.e. 


^2/^  7^P2/^2  =  9i/2  +  — ^Pl/Pl 


^~^2i  2  ^~^2i  2 

——qj  +  fla  =  — ;r-9i  +  “i 


-  2  * 


(3) 


where  {J^j}  denotes  sonic  condition  in  the  upstream  flow,  when  achieved  isen- 
tropically. 


The  set  of  four  equations  ( 1)  -  (3)  yields  the  shock  wave  exit  conditions,  qa,  p2^ 
Pa  and  aa,  in  terms  of  the  corresponding  entry  conditions  and  of  q^.  The 
solution  is  obtained  as  follows: 


Dividing  (2)  by  (1)  and  substituting  for  and  af  from  (3)  yields 


92  +  q^/q2  =  ?i  +  9»/9i 


or  using  (1)  again 

92  +  (9«/^i9i)P2  =  91  + 

Adopting  *  to  denote  that  the  speeds  and  densities  in  (1)  and  (4)  have  been 
normalized  by  divisors  q^  and  and  dropping  the  suffix  “2”,  gives 


p4  =  Pi4i 


and 


(!') 


9  +  ^/(Pi9i)  =  91  +  1/91 


(4') 
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Thus  tor  a  given  shock  wave  entry  condition  (qi./Ji),  referring  to  Fig.  .5  the 
shock  wave  exit  condition  (q2,P2)  <s  represented  by  the  subsonic  intersection  B 
of  rectangular  hyperbola  ( T)  (the  continuity  condition)  with  straight  line  [A') 
in  the  (q.p)  plane,  the  other  (supersonic)  intersection  being  obviously  the  shock 
wave  entry  point  A  (qi.pi)  in  Fig.  .5. 

Eliminating  p  from  equations  {!')  and  (4')  yields  a  quadratic  in  q; — 

-  (gi  +  1/91)9  +1  =  0 

with  roots  q2  and  qi.  Hence  the  product  of  the  roots  satisfies 

9291  =  1  l'^) 

and  combining  this  with  (1)  normalized  gives 

h/P\  =  9i  (6) 

With  the  aid  of  equations  (5)  and  (6)  it  is  possible  to  derive  an  analytical  p  vs. 
q  relationship  at  the  exit  from  a  normal  shock  wave.  i.e.  to  derive  the  locus  of 
point  “B”  in  Fig.  5  as  the  parameters  qt  and  pi  of  the  shock  entry  condition 
vary.  This  is  achieved  as  follows: 

.At  shock  entry  (suffix  “1”)  the  flow,  which  is  isenergic  throughout,  has  up  to 
that  point  been  also  isentropic.  So  substituting  pi  and  qi  for  p  and  q  in  the 
isentropic,  isenergic  channel  flow  equation 


leads  to 


With  (6)  this  gives 


7  -  1  ,2^ 
2 


h  =  9? 


/7  +  I  7-1  .2\1^ 

(2  2  ’V 


and  using  (5) 


i.e.  dropping  suffix  “2” 


(•) 

(8) 


(9) 


« 
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whicn  simplifies  to 


where 


R  =  \/e 


O') 


The  shock  exit  (density  vs.  speed)  relationship  (9)  is  presented  on  Fig.  5  the 
Pex  '"S-  q  curve,  for  l>q>  As  the  normalized  shock  e.xit  speed  q 

decreases  from  1.0,  p^x  at  first  increaises  up  to  a  maximum  value  of  about  1.17; 
it  then  falls  sharply  to  zero  at  q=^- 

The  maximum  occurs  when  ^  =  0  that  is.  differentiating  (9')  logarithmi¬ 
cally  w.r.t.  R,  when 


l/R- 


=  0 


giving  the  normalized  exit  speed 


q  =  1/yS  = 


.7638 


(corresponding  to  a  shock  inlet  Mach  number  of  y/2).  Substituting  this  q  value 
into  (9)  yields 

{pex)MAX  =  I -16605 

Shown  also  on  Fig.  5  is  the  isentropic,  isenergic  flow  relationship  (7)  presented 
as  the  Pit  vs.  q  curve,  for  the  complete  q  range  from  (=  VS  and  corre¬ 

sponds  to  infinite  Mach  number)  down  to  0.  The  subsonic  intersection  C  of 
this  curve  with  the  continuity  equation  (1')  curve  gives  the  shock  exit  conditions 
as  predicted  by  the  conventional  FPE  method.  The  wide  gap  between  inter¬ 
sections  B  and  C  thus  represents  the  difference  between  the  Rankine-Hugoniot 
and  FPE  shock  exit  conditions.  It  demonstrates  again,  quite  dramatically, 
the  large  and  rapidly  increasing  error  {pc  -  pb)  in  the  shock  exit  density  when 
determined  by  the  conventional  FPE  method,  particularly  when  B  is  to  the  left 
of  peak  “n”  of  the  /3«x  curve  (corresponding  to  the  shock  inlet  Mach  number 
exceeding  v'^)- 

So  within  a  normal  shock  wave  the  following  two-part  (density~velocity)  rela¬ 
tionship  will  ensure  Rankine-Hugoniot  conditions  at  shock  exit: — 


q  >  1  p  =  pit(q)  -  equation  (7) 

q  <  I  p  =  Pex{q)  -  equation  (9) 


Adopting  this  approach,  the  only  modification  required  to  the  FPE  flow  model 
within  the  shock  wave  is  to  use  equation  (9)  instead  of  equation  (7)  for  the 
subsonic  (q<  1)  portion  of  the  wave. 


* 
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Note  that  equation  (9)  is  conveniently  parameter-free.  i.e.  it  requires  no  knowl¬ 
edge  of  the  shock  wave  inlet  conditions  qj  and  pi  (unlike  equation  (4')). 

Flow  conditions  downstream  of  the  shock,  consistent  with  the  modified  shock 
wave  model,  are  dealt  with  later  at  §4. 

The  use  of  equation  (9)  has,  however,  two  potential  drawbacks,  cdthough  these 
may  be  only  slight  or  non-existent  in  many  cases.  Firstly  the  exit  stagnation 
density,  which  is  pivotal  in  determining  the  flow  downstream  of  the  shock  wave, 
is  affected  by  any  error  in  the  perceived  shock  exit  location  and  state.  Secondly, 
as  indicated  above  and  in  Fig.  5,  equation  (9J  when  plotted  hais  a  slight  peak 
—  i.e.  ceases  to  be  monotonic  —  if  q  <  .7638;  this  happens  only  if  Min  >  v/2. 
The  peak,  if  it  occurs,  does  so  towards  the  shock  wave  exit  and  the  density  then 
declines  until  shock  exit  is  reached.  For  Min  <  1-6  the  decline  is  less  than  2t%. 
Some  local  numerical  instability  might  conceivably  result  from  this  peaking, 
during  a  CFD  iterative  solution  of  a  transonic  flow.  However  such  instability 
could  be  prevented  by  clipping  the  density  peak,  so  replacing  it  with  a  density 
plateau,  and  then  progressively  eliminating  the  clipping  in  suitably  conrolled 
stages  as  the  iterative  solution  proceeds. 

To  avoid  these  potential  problems  another  (density~velocity)  shock  wave  wafer 
relationship  is  developed  in  the  following  section. 


3  SECOND  PROPOSED  DENSITY  vs.  VELOC¬ 
ITY  RELATIONSHIP  WITHIN  THE  SHOCK 
WAVE  WAFER 


Referring  again  to  Fig.  5,  virtually  any  more-or-less  monotonic  curve  joining 
A  to  B  (other  than  the  continuity  equation  {!')  curve)  appears  suitable  as 
a  replacement  here  for  the  isentropic,  isenergic  flow  equation  (7)  of  the  FPE 
method.  Such  a  curve  would  intersect  the  continuity  equation  ( 1')  subsonically 
at  B,  which  is  on  the  (p„  vs.  q)  curve  —  thereby  ensuring  Rankine-Hugoniot 
exit  conditions. 

However,  as  previously  mentioned,  there  exists  a  possible  problem  concerned 
with  pinpointing  computationally  the  precise  location  of  X,  the  exit  point  of 
a  streamline  from  the  shock  wave  wafer  (corresponding  to  state  B  in  Fig.  5).* 
This  is  another  changeover  point  for  the  (density  vs.  speed)  relationship:  pro¬ 
ceeding  streamwise  along  the  streamline  through  X  the  (p  vs.  q)  relationship 
changes  there  from  the  shocx  wave  wafer  internal  form 

^One  way  of  estabiishing  computationally  the  location  of  X  is  to  compute  the  spatial  density 
gradient  with  respect  to  the  streamsrise,  or  nearly  streamwise.  coordinate  of  the  flowfield.  In 
the  immediate  vicinity  of  X  this  gradient  is  extremely  large  just  upstream  of  X  (being  within 
the  shock  wave  wafer)  but  very  much  smaller  just  downstream  of  X.  Such  a  test  thus  enables 
.X  to  be  located  approximately,  with  fair  accuracy. 
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(not  yet  determined)  to  the  isentropic,^  isenergic  form  for  flow  downstream  of 


Although  this  change  is  piece-wise  continuous  in  ordinate,  it  is  not  necessarily 
so  in  gradient.  If  due  to  the  computational  discretization  of  the  flow  field, 
or  for  any  other  reason,  a  small  normalized  speed  error  Aq  arises  at  X  from 
a  slightly  inaccurate  determination  of  the  precise  shock  exit  location  X.  then 


there  will  also  be  an  error 


in  the  value  of  pg  ( used  in  (11)).  Consequently  at  the  presumed  (but  slightly 
incorrect)  shock  exit  state  B  there  will  be  a  fractional  error  —  either  overshoot 
or  undershoot  —  in  the  normalized  stagnation  density 


equal  to 


«  -  (  ^  ^Psiu  ,  _ 24 _ N 

dq  7  + I -(7-  1)9Vb 


(12) 


correct  to  the  first  order  in  Aq.  As  an  accurate  computation  of  the  flow-field 
downstream  of  the  shock  wave  depends  on  accurately  determining  the  stagna¬ 
tion  density  at  shock  wave  exit^,  it  is  important  that  the  form  of  psw(d)  should 
be  chosen  to  make  this  error  zero,  if  possible,  but  otherwise  minimal. 

Therefore  from  equation  (12)  a  further  desired  condition,  additional  to  those 
listed  at  §1,  is 

fj_d^\  ^  f - 24 - 

\Psw  dq  Jg  \7  +  1  -  (7  -  1  )9Vs 
Equation  (13)  in  fact  expresses  the  criterion  that  the  stagnation  density  —  as 
derived  from  (10)  —  is  to  be  stationary  just  upstream  of  X.  Downstream  from 
X  the  stagnation  density  of  the  flow  is  constant  along  a  given  streamline  (as 
expressed  at  equation  (11)).  Thus  along  a  streamline,  by  imposing  condition 
(13)  there  will  be  at  X  piece- wise  continuity  of  the  stagnation  density  not  just 
in  ordinate  but  also  in  first  derivative  w.r.t.  q. 

This  piece-wise  continuity  of  stagnation  density  at  X  could  be  extended  to  any 
desired  order  of  derivative  w.r.t.  q  by  suitable  choice  of  the  function  Psw(q)  'u 
( 10).  However  it  is  best  to  restrict  (10)  to  a  second  order  polynomial  in  q;  third 
and  higher  order  polynomials  run  the  risk  of  possible  non-physical  intersections 
—  i.e.  intersections  other  than  A  and  B  in  Fig.  5  —  with  equation  (1'). 


^.More  precisely,  isentropic  along  a  streamline. 
'Refer  to  equation  (11). 


From  (  5 )  and  (6)  B  (qj,  p^)  is  ( l/qj.  Piqi^);  A  is  (qj,  pi ).  Also,  substituting 
these  B  coordinates  into  (13),  the  desired  slope  at  B  of  the  proposed 
quadratic  relationship  (101  is 

-2pw?/[(7+  1)9?  -(7  -  1)1 


The  desired  (p,q)  equation  (10)  within  the  shock  wave  wafer  is  then  readily 
verified  as  being 


P/Pi  =  (I  +9?)  -  9‘i9  -  (9-9i)(9-  l/9i)9?/ 


7-1 


(14) 


Note  that  equation  (14)  really  contains  only  one  parameter,  qi,  because  px 
is  obtained  immediately  in  terms  of  qi  from  equation  (8).  Knowledge  of  qj 
requires  the  determination  of  shock  wave  entry  location;  this  is  achieved  in 
similar  fashion  to  determining  the  shock  wave  exit  X.  (See  footnote  2.) 

An  example  of  equation  (14),  for  some  chosen  value  of  qi,  is  plotted  in  Fig.  6. 


3.1  INTERSECTIONS  OF  EQUATION  (14)  WITH  THE  CON- 
TINUITY  EQUATION 


Equation(14)  intersects  the  continuity  equation  (1')  at  just  three  points:  for 
substituting  (T)  into  (14)  yields  a  cubic  in  q: — 

4i  =  (l  +  9x)9-9i9*-9(9-9i)(9-l/9i)Qi  (15) 


where 


Qi 


9?/ 


7+1 

7-1 


Note  here  that  since  qi  >  1,  therefore  Qi  >  ^  >  0. 
( 15)  is: — 


(16) 


-Qig^  +  [Qi(9i  +  l/9i)  -  9i]9^  +  (9?  +  1  “  Qi)9  -  9i  -  0 


which  factorizes  to 


(-<3i9  -9i)(9-  9i)(9-  1/9i)  =  0 


The  roots  of  this  are 

9  =  9i  (point  A), 

9  =  1/91  (point  B) 

and  9  =  -qi/Qi 

Now  Qi  >  0.  (Refer  below  (16).)  Hence  the  third  intersection  (q=  -qi/Qi) 
is  illusory  because  it  corresponds  to  a  negative  value  of  q.  Thus  there  are  only- 
two  possible  intersections  of  (14)  and  (T)  —  the  two  physical  ones,  .\  and  B. 
as  required.  (See  Fig.  6.) 
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3.2  MONOTONICITY  OF  EQUATION  (14) 


It  is  readily  verified  that  equation  (14)  may  be  written  cis 


P/Pl  = 

+  [(Qi  -  l)4i  +  Qt/?il4  +  [4i  +  1  ~  Qi] 

(14' 

Differentiating, 

^  ^  =  -2Q  i4  +  [((?  1  -  1)91  +  Q 1  /9i  ] 

and 

d^d  - 

—  =  -IQiPi  <  0 


So  p  is  stationary  (a  maximum  in  fact)  only  at 

9  =  4  =  [(1- i/«i)4i  +  i/4i]/2  =  i/4i-[i/4i-(i~i/«i)4i]/2 

=  from  (16) 


Thus,  qi  being  positive,  q  <l/qi  (=  qj,  the  q  value  at  B).  Hence  in  Fig.  6 
the  vertex  of  parabola  (14')  is  to  the  left  of  B.  Therefore  over  the  range  of 
application  A  to  B  parabolic  arc  (14')  is  monotonic  and  convex  upward;  within 
this  range  (as  shown  previously)  it  intersects  the  continuity  equation  ( 1')  only 
at  A  and  B. 


4  CONDITIONS  DOWNSTREAM  OF  THE  SHOCK 
WAVE 

In  the  FPE  flow  field  model  the  flow  is  assumed  to  be  universally  isentropic; 
therefore  the  stagnation  density  and  stagnation  pressure  are  globally  constant 
throughout  the  flow  field.  With  the  modified  flow  field  model  resulting  from  the 
new  shock  wave  representation,  however,  these  quantities  —  being  dependent 
on  the  local  shock  inlet  Mach  number  Mi  —  now  vary  along  the  exit  face  of  an 
embedded  normal  shock  wave  from  its  foot  ( Mi  =body  surface  Mi  >  1)  to  its 
tip  (Mi  =  l),  and  then  remain  constant  downstream  of  the  shock  only  along  a 
given  streamline. 


4.1  DETERMINING  THE  STAGNATION  DENSITY  AT  SHOCK 
WAVE  EXIT 


Stagnation  density  is  defined  by 

PSTAG  =  p! 


7- 

1 - -q*' 

7  +  1 


Therefore  at  shock  wave  exit  the  normalized  stagnation  density  PstaGb  1®  given 
bv 


hjAGg  =  PSTAGbIpIo  =  pb/ 


1  '1'  "  ^  -2 
‘■TTI’® 


1/h-i) 
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However  exit  state  B  lies  on  the  B-locus  (p  vs.  q)  curve  defined  by  equation 

(9).  So 


PB  = 


'7  +  1-2 


i/h-i) 


Therefore  substituting  for  pB, 


PSTAGb  = 


(!■) 


4.2  BODY  SURFACE  PRESSURES 


The  normalized  stagnation  pressure  at  shock  wave  exit,  psTAGgi  readily 
derived  from  the  fact  that  the  stagnation  temperature,  Tstag.  is  constant 
throughout  the  isenergic  flow.  Thus: — 


PSTAGb  _  PSTAGb  '^STAGb  _  PSTAGb 
PSTAGoa  PSTAGcc  TsTAGoo  PSTAG^ 


(18) 


(where  pressures  are  normalized  by  the  divisor  pJo). 
Now 


PSTAGb  =  PstaGooIpIo  =  critical  density  ratio  = 


l/(-r-l) 


and 


PSTAGb  =  PstaGooIpIo  =  critical  presssure  ratio  =  ) 


V(-r-l) 


Therefore  inserting  these  values  into  (18) 


PSTAGb 


PSTAGb 


with  PSTAGb  determined  from  (17),  Along  the  body  surface  streamline  pstag 
is  (discontinuously)  constant,  with  the  values  pstaGoo  upstream  and  pstaGb 
downstream,  respectively,  of  the  shock  wave.  Also  q  is  known  throughout  the 
flowfield,  and  in  particular  on  the  body  surface,  from  the  converged  flow  solu¬ 
tion.  Hence  the  body  surface  pressure  distribution  is  immediately  calculable 
from 


P  =  PSTAG 


I 


7  -  1  -2 

7+  l’ 


7/(7- 1 ) 


5  CONCLUDING  REMARKS 

1.  In  view  of  the  importance  of  correctly  predicting  shock  wave  exit  conditions 
in  CFD  solutions  of  shocked  transonic  flows,  a  new  improved  relationship  has 


* 
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been  'lerived  between  density  (p)  and  velocity  (q)  inside  an  inviscid  normal 
shock  wave.  This  relationship  is  much  superior  to  that  conventionally  applied 
there  by  the  conservative  full  potential  equation  (FPE)  method  and  should 
improve  very  considerably  its  accuracy  and  range  of  applicability,  especially  for 
flows  containing  strong  shock  waves. 

2.  In  deriving  the  improved  (p,q)  relationship  the  FPE  assumption  of  globally 
isentropic  flow  has  been  discarded  inside  the  normal  shock  wave  and.  instead, 
the  correct  physical  conditions  at  shock  exit  have  been  enforced.  .Additionally, 
the  relationship  (equation  (14))  ensures  that  the  shock  exit  stagnation  density, 
a  vital  element  in  determining  the  post-sho«dt  flow,  is  unaffected  (to  the  first 
order)  by  any  small  error  that  may  be  present  in  the  perceived  shock  wave  e.xit 
location  (and  state). 

3.  The  only  parameter  appearing  in  the  new  (p,q)  relationship  (14)  is  the 
normalized  velocity  (or  alternatively  the  Mach  number)  at  shock  entry.  This 
can  be  determined  —  and  continuously  updated  —  from  the  developing  iterative 
flow  field  solution  once  the  shock  becomes  recognisable. 

4.  An  alternative  version  of  the  (p,q)  relationship,  equation  (9),  which  has  the 
advantage  of  being  parameter-free,  is  also  presented.  With  this  relationship, 
however,  the  exit  stagnation  density  is  subject  to  any  error  in  the  perceived 
shock  wave  exit  location  (and  state).  As  weU,  if  the  normal  shock  wave  inlet 
Mach  number  Mj  is  higher  than  y/2  a  slight  non-monotonicity  develops  in  this 
alternative  (p,q)  relationship.  This  occurs  near  to  the  shock  wave  exit,  with 
the  density  peaking  and  then  declining  slightly  to  its  value  at  exit.  The  den¬ 
sity  decline  is  at  most  only  2^%  for  Mi  <1.6.  In  theory  this  peaking  effect 
might  conceivably  trigger  some  numerical  instability  in  an  iterative  computa¬ 
tional  solution.  Such  instability  can  be  avoided  by  clipping  the  density  peak, 
thus  replacing  it  with  a  density  plateau,  and  then  progressively  eliminating  the 
clipping  in  suitably  controlled  stages  as  the  iterative  solution  proceeds. 

.5.  The  two  versions  of  the  improved  (p,q)  relationship  applicable  inside  an 
inviscid  normal  shock  wave  are  given  by  equations  (14)  and  (9)  (as  stated); 
CFD  usage  will  indicate  which  is  superior. 

6.  .A  simple  rule  for  extending  the  (p.q)  relationship  to  the  case  of  an  oblique 
shock  wave  is  given  in  the  .Appendix:  equations  (.A22),  (.A22ob),  (.A25)  and 
(.A26)  refer. 
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Appendix  A 


EXTENSION  TO  OBLIQUE  SHOCK  WAVES 


Figure  Al:  OBLIQUE  SHOCKWAVE  STATES 
A1  FUNDAMENTAL  AND  DERIVED  EQUATIONS 


When  the  shock  wave  is  oblique,  the  continuity,  momentum  and  energy  equa¬ 
tions  become  with  the  notation  of  Fig.  Al: — 

continuity 

P29jn  =  Piqin  =  mflux  (Al) 


(where  mflux  is  the  mass  flux  rate  normal  to  the  shock  wave  and  (qin.qt)> 
(qjntqt)  are  the  velocity  components  normal  and  parallel  to  the  shockwave  at 
entry  and  exit) 


momentum 

{P2 -Pi)  =  mflux(gj„ -<7jn)  (A2) 


k 
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which,  using  equation  may  be  written' as 

VPi  PiJ 


2  2  _  (<7ln  <72n)(P2  -  Pi  i 

9ln  “  *72n  -  - 


mflux 


(P2  -  Pi  1 


energy  (  §  'I') 


energy  (  §  '2' ) 


SLl  4.Sl]  ^  SL-El  =  >  +  ^  .2 

•2  -2  7-lpi  ■2(7-1)’’° 


P2 


T'  +  1  .2 

?oo 


2  2;  7-lP2  -2(7-1) 


I  A3 ) 


(A4) 


(Ao) 


Expressions  for  qjn  and  p-i  analagous  to  those  for  the  normal  shock  wave  case 
are  found  as  follows. 


Differencing  equations  (A4)  and  (A5)  gives 

2  2  27  /  P2  Pi  \ 

’--’2n  =  —  (---j 

Combining  equations  (A3)  and  (A6) 

\Px  P^)^  7-iVp2  Pi/ 


(A6) 


i.e. 


which  leads  to 


\Pl  /  \Pl  /  7-1  Vpi  PiJ 

( 7±i')  £1  _  I 

P2  V-»-i/ei  '■ 


Pi 

\y-l  /  Pi 

Pi/ Pi  may  be  found  in  terms  of  the  initial  state  '1'  as  follows;- 
From  equations  (  A2)  and  (Al) 

Pi  Pi 


(A7) 


Qln  -  qin  = 


Piqin  Piqin 


=  (1^) 


-1-1-1  -2  ll 

''-*•1  .  ?i. 

] 

[2(->-l)^^  2 

qin 

2(-,-liPoo  2 

<  '  —  ■ 

^2n 

1 

2 

Pin 

M 

from  equations  (A5)  and  (A4).  Hence 


(Pin  -  qin)  = 


(Pin  -  qin) 


\r-^)q'J  -  i^)qf 

qinqin 


+ 


y  —  1 


Since  q2n  ^qin-  therefore 


:^)q^-(^)<if  ^7-1 

P2nPln  7 


* 
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giving 


From  (Al)  and  (A8) 


92n9ln  ■ 

=  <?»-! 

1  A^ 

Pi  _  9ln 

9ln 

cr> 

< 

Pi  92n 

9*  - 

So  with  equations  (AS)  and  {A9)  —  together  with  energy  equations  i  A4)  and 
( A5)  -  in  mind,  replace  q;^  as  a  normalizer  by  q*  =  -  t  ^)qtV  Equation 

(.A8)  then  becomes 

which  is  analogous  to  equation  (5)  of  the  main  text  for  a  normal  shock  wave. 
Now  referring  to  rig.  Al, 

=  [9:0*  -  (^)9?C0S^^]  (All) 

3  being  the  acute  angle  between  the  obUque  shock  wave  and  the  velocity  vector 
qi  at  entry 

’  (^)  [(’  ^ 


=  9« 


(from  main  text  equation 
(3)  with  Mi=qi/ai) 


[  ^  (-r-nM? 

(using  (3)  again). 


i.e.  since  qisind  =  qj^ 

^.2  _  f,  ,  2  .  Jl  +  ( ^).V/2sir 


A 12) 


A2  COMPARISON  WITH  NORMAL  SHOCK  WAVE  RE¬ 
SULTS 


.Vow  from  (.All)  the  energy  equations  (A4)  and  (.A5)  become:- 

SL  +  _2_£l  =  JL±1_5-2 

2  7-lpi  2(7-1)’ 


and 


Sh  +  —l_El  =  JL±_L-2 

2  7-lP2  2(7-1)’ 


(.A4') 


(.\5') 
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Also  the  continuity  and  momentum  equations  ( previousiy  presented  in  this  Ap¬ 
pendix!  are;  — 

Pl9ln  =  /»292n  =  mfluX  (Al) 


and 


iP7-  P\)  =  mflux(gi„  -  q-2n) 


(A2) 


These  four  equations  are  formally  identical  to  those  for  a  normal  shock  wave 
provided  that: — 


Eqn.  nos. 
(A1),(A2),(A4^ 

(Al),  (A2),  (A5') 

(A12) 


OBLIQUE  S.W. 
gin 


q2n 


— >  Replacing  — 

H 


It 


n 


NORMAL  S.W. 

qi 

q2 


From  the  above  four  basic  equations  (A4'),  (A5').  (Al)  and  (.A2)  the  following 
results  can  be  derived  (via  the  equations  tabulated): — 


(AlO) 

qinq2ii=q*^ 

// 

qiq2=q® 

(A9)  with  (A12) 

PI  ~  2+(-y-l)MJ»u»^p 

n 

p,  (-r+DMf 

PI  ~  2+(-r-l)MJ 

DITTO  with  (A7) 

n 

E2  — 

^  -  1+1 

{qin  and  q2n  replace  qi  and  q2 
q’  replaces  q^ 

Misin/d  replaces  Mi 


(A13) 


then  the  basic  equations  and  the  results  derived  therefrom  are  the  same  for  an 
oblique  shock  wave  as  for  a  normal  shock  wave. 


One  further  replacement  is  required,  that  for  the  density  normalizer  This 
is  obtained  as  follows: — 


Considering  the  density  at  the  shock  wave  entry  face 


7  -  1 
7  +  1 


7  -  1 
7  +  1 


(A14) 


(This  is  essentially  the  same  equation  as  m«un  text  equation  (8).)  Therefore 


2 


=  7  +  1  -  (7  -  1) 


(A15) 


<71  =  <7in  +  <7?  \ 

+  J 


Now 

and  from  (.All) 
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(A16) 


Substituting  from  ( A16)  for  qf  and  into  (Alo)  gives 


=  [(7  +  1)?*^  -  (7  -  1  )«??„]  /<? 


Therefore 


pt  ‘ ^  '  '’.’./’-I 


Now  let 


Therefore  from  (A  12) 


P*  =  P»(9*/<7»)'^‘^-‘’ 


(A18) 


j’ =  P«{[l  +  ^-WiW3  /  l  +  (A  19) 


Substituting  (A18)  into  (A17)  leads  to 

(ei)’"  =  f, .  /  fi .  l::!] 

\P-)  [  1  +  l\?Vj'l  7  +  lJ 

with  p*  given  by  (A18). 

(A20)  is  formadly  identical  to  main  text  equation  (8)  but  with 

!p’  replacing  \ 

q*  replacing  q^^  >  1 

qi„  replacing  qi  ) 

Thus  for  an  oblique  shock  wave  the  appropriate  density  normalizer  is  p'. 


(A20) 


(A21) 


A3  DENSITY  vs.  VELOCITY  RELATIONSHIP  FOR  OBLIQUE 
SHOCKS 


The  ip  vs.  q)  law  that  has  been  derived  for  the  interior  of  a  normal  shock  wave 
can  be  written  as  ,  . 

P  f. _ _ I  9  1  ! 


=  function 


ion  1  — r 

W) 


(A22) 


This  equation  is  adapted  to  the  oblique  shock  wave  case  by  applying  appropriate 
replacements  as  indicated  by  (A13)  and  (A21),  yielding 


^  =  function  (-%)  =  function 

P' 


ion  — 7  -  T-:r  = 

q-i) 


function! -  Aja) 
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p  =  A;tfunction(fc2^^  -  ^3) 


( A'2'2ob ) 


i.e. 

with: 


=  P’Ip\o  =  {[l+0.2.W,W;3)/[l  +  0.2.U?]}2-5 

^2  =  =  [1 +  0.2Afi*l/[l  +  0.2A/2sin2j] 


*=3  =  =  {qiC0&^l3/q'^)(q^/q’^) 


■^A/?cos2/3‘ 

.  1  +  ^A/2_ 

1  + 

by  main  text  equation  (3)  and  equation  (A12); 


hence 


by  (A  19) 
by  (A  12) 


fca  =  1.2A/j*cos^/3/[l  +  0.2A/jsin*/J]  J 

(A23) 

Using  the  isentropic  isenergic  flow  relationship  at  equation  (3)  again.  .Mj  is 
given  by 

-^i  =  5(9?/9«')/{6  -  =  59?/(6  -  9?)  ( A24) 


leading  to 

.  ,  i 

*  1  -  (qJcos*/3)/6 

From  (A23)  ki,  kj  and  k3  are  interrelated  thus: 

*1  =  1 

and  kj  =  6{*2  ~  1)  / 


(A25) 


(A26) 


Equation  (A22)  thus  converts,  by  simple  linear  transformations,  into  the  (more 
general)  oblique  shock  wave  case  (A22ob),  with  the  aid  of  equations  (.A25)  and 
(.A26). 


Note  that  the  angle  f3  appearing  in  equation  (.A25)  is  a  function  of  Mi  and  6. 
the  flow  deflection  angle  (see  Fig.  Al).  Hence  from  (,A25)  and  (.A24)  k2  is  a 
function  of  qi  and  6. 
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Fig  3.  Pressure  Ratios  for  Normal  Shocks 


Fig  5.  {p  vs  q)  Relationships  for  Isentropic  Flow  and  for  Normal  Shock  Waves 
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